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Abstract 

We derive an explicit expression for the quantum conductivity of 
a molecular wire containing a redox center, which is embedded in 
an electrochemical environment. The redox center interacts with the 
solvent, and the average over the solvent configurations is performed 
numerically. Explicit calculations have been performed for a chain of 
three atoms. When the redox center interacts strongly with neigh- 
boring electronic levels, the current-potential curves show interesting 
features like rectification, current plateaus and negative differential 
resistance. Electronic spectroscopy of intermediate states can be per- 
formed at constant small bias by varying the electrochemical potential 
of the wire. 



1 Introduction 

Understanding of electron transport through a single molecules received an 
increased interest due to the speculation of employing molecular units as fun- 
damental elements of computer circuits [1,2]. Additionally, electron transfer 
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in molecular wires at nanoscale level received further attention, both at the 
level of formalism as well as ab-initio calculations, due to its possible relevence 
in understanding and application for a class of diverse problems like sensors, 
photonics, solar energy conversion [3] . Controlled charge movement in a suit- 
ably designed molecule can be used as the basis for storing and processing 
of information. Such Quantum-Dot Cellular Automata architecture has been 
experimentally realised in a series of experiments for a variety of applicable 
components like memory cells, logic gates and clocked memory cells [4-6]. 
A practical implementation of QCA architecture consists of a single redox 
center with an organic or inorganic bridging group [7] . The efficiency of solar 
energy conversion process depends not only on efficient photon capture but 
also on charge seperation and transport through very large distances. Since 
the charges are created by sunlight on the surface of an assembly of molecules 
or semiconductors, it is resonable to expect molecular wires to act as rele- 
vant acceptors of the charges. Also the weak solar fluxes imply a very low 
current and a need for fast charge transport. Certain classes of polymers and 
oligomers have been proposed as ideal candidates for satisfying the above 
criterion for increasing the yield in solar energy conversion [8,9]. 

The above mentioned are some of the reasons for the recent increased 
surge in interest for understanding charge transport along molecular wires. 
Typical theoretical work in the field involves obtaining generic expressions 
for the conductance, current- volt age profiles, rate constants, transition proba- 
bilites etc. Formal works on transport properties along molecular wires were 
carried out on Donor-Bridge-Acceptor (DBA complexes) systems, wherein 
electrons are transferred between donor and acceptor connected by a molec- 
ular bridge [10-13]. Further, electron transfer between reservoirs connected 
by a molecular bridge has been studied by Ratner and co-workers [14-16]. 
The above works resort to time- dependent quantum mechanics for obtain- 
ing expressions relevant to electron transfer. Expressions for the conductance 
between two reservoirs connected by monoatomic sites were well known in 
mesoscopic physics. The formal expression was first derived by Caroli and 
co-workers [17] and was later expanded to a broader class of problems by 
Wingreen and co-workers [18-20]. Recently, the same expressions were re- 
derived by various authors [21-25] by formulating Quantum Langevin Equa- 
tions (QLE). Initially, both the molecular wire and mesoscopic conduction 
were modelled using tight-binding Hamiltonians, and since at a Hamilto- 
nian level these problems seem identical, it is expected that the expression 
obtained for one should be applicable for the other. 
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Several authors have pursued other computational methods such as den- 
sity functional theory, first principle ab-initio calculations and package simu- 
lation of Non-Equilibrium Green's Function in studying conduction through 
molecular wires [26-30]. Most of these authors differ in their treatment of 
the electrodes and the interaction of the metal-molecule couphng. Most of 
the earlier works were oriented towards a better approximation for modelling 
the interaction and self-energies at the molecule-metal junctions, while little 
work has been done in incorporating the effect of additional interactions the 
electron might experience in the molecular wires. Though numerous works 
have been done of the subject of electron-phonon coupling with relevance to 
quantum dots [18,31-36], including treatments for classical, quantum, equi- 
libriated and out of equilibrium systems, exact treatment of such a process 
in specific for molecular wire has not received much attention. Molecular 
wires differ from quantum dots in that the observed conductance behaviour 
of quantum dot is dominated by Coloumb blockade. 

The present work focuses on a special electrochemical case: a molecular 
wire containing a redox-center connecting two electrodes. The electrochemi- 
cal case is of special interest since two potentials can be varied independently: 
the bias between the two electrodes, and the potential of one of the elec- 
trodes with respect to the solution. The latter acts hke a gate voltage that 
controls the current in the wire. In addition, the redox center interacts with 
the solvent, whose fluctuation will affect the current. The specialties of the 
electrochemical situation were first elucidated by theorists [37-40]. Starting 
with the pioneering paper of Tao [41], there have been a fair number of exper- 
imental studies of electrochemical systems [42-48,50,51] which in turn have 
generated more theoretical work (see e.g. [52,53] and references therein). 

Most of the theoretical work on electrochemical systems has been re- 
stricted to special systems with one or two intervening redox centers. In this 
study, we will consider a wire of arbitrary lengths containing one redox center 
interacting with the solvent. Using a tight-binding Hamiltonian and Green's 
function techniques we will derive an expression for the current which is ex- 
act for the case where the interaction to the two electrodes can be treated in 
the wide-band approximation. These calculations will be illustrated by model 
calculations for particularly interesting cases: steps and negative differential 
resistance, and spectroscopy of intermediate electronic states. 
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2 Model Hamiltonian, Green's functions and 
current density 

The model system that we consider consists of two metal electrodes, labeled 
R and L, connected by a chain of 2n + 1 atoms - an odd number is chosen 
for convenience only. The atom in the center is redox-active and interacts 
with the solvent; thus, we identify the index n + 1 with the index r of the 
redox species. We use a tight-binding model, in which each atom contains 
one orbital and interacts only with its nearest neighbor. The corresponding 
Hamiltonian can be written in the form: 



2n+l 



H = ^ ^ ek,ink,i + trTLr + e 

k i=L,R i=l 
n— 1 2n 

+ ^{vic\ci+i + h.c} + ^ {vic\ci+i + h.c} 

i=l i=n+2 

+ ^ {Vic\ci+i + h.c} + ^{^fc,l4,LCl + ^fc,2n+l4,RC2„+l + h.c} 

i=n,n+l{r) k 

+ ^ ^ hu]^ql + ^ hw^g^q^Ur (1) 

In this Hamiltonian, n always denotes an occupation number, a cre- 
ation c an annihilation operator, e an energy, and v a coupling constant. 
The first line contains the diagonal elements, the indices (/c, L) and (fc, R) 
labeling the electronic states on the two electrodes. The second and third 
lines give hopping elements between adjacent sites, and the last line the po- 
tential energy of the solvent, with coordinates q^, and frequencies and its 
interaction with the redox center r; the are the corresponding coupling 
constants. Equation ([1]) is a natural generalization of the Hamiltonian for 
redox-mediated tunneling via one center [40] . 

The matrix form of the fermionic part of the above Hamiltonian Hp has 
the generic form: 
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Working within the tight-binding model, it is understood that Hp has 
non-zero entries only in diagonal and sub-diagonal elements. The general 
scheme of approach is to calculate the quantum conductance and then obtain 
the current by integrating the conductance between appropriate limits. The 
formula employed for obtaining the quantum conductance, or tunneling rate, 
is the same as the one used by Datta et al. [55]. This form of the formula 
was first derived by Caroli et al. [17] and was subsequently derived in a much 
wider context by Wingreen and co-workers [18-20] . 



g ^ Tr[G^riG«rR] (2) 

As before, the subscripts L and R refer to the left and right reservoir. F 
denotes the imaginary part of self-energy (for ease of following the notations, 
r = Vk{ImG\j^)vl). The quantity of interest \s \ {I \ G \ 2n-\- 1) |^, where G 
is the Green's function obtained from the above Hamiltonian. This can be 
obtained by separating the Hamiltonian into two parts: H = Hq + V and 
considering a Dyson equation , G = G^ + G^VG, where the simplification of 
the problem results from the choice of V. Letting V — Yl:i=n n+i ^ic\cij^i + h.c, 
the Hamiltonian contains 3 block matrices. Physically this amounts to 
cutting the 2n + 1 atom chain at 2 places on either side of the redox couple. 
The closed form for the element (1 | G | 2n -|- 1) is obtained as shown: 

(1 I G I 2n + l) = (1 I G° I 2n + l) + ^(l | G° | \ V \ \G\2n + l) 

(3) 

= (1 I I n)(n I 1/ I n + l)(n + 1 I G I 2n + 1) 
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(n + 1 I G I 2n + 1) = (n + 1 I I n + l)(n + 1 | 1/ | n)(n | G | 2n + 1) 

(4) 

+ (n + 1 I G° I n+ l)(n + 1 | V | n + 2)(n + 2 | G | 2n + 1) 



(n I G I 2n + l) = (n I G° I n){n \ V | n+ + l | G | 2n + 1) (5) 



(n + 2 I G I 2n + 1) = (n + 2 I G° I 2n + 1) (6) 
+ (n + 2 I G° I n + 2)(n + 2 | V | n+ l)(n+ 1 | G | 2n + 1) 

Prom the above 4 equations , Gi,2n+i = (1 | G | 2n + 1) can be solved: 

1 — G^^j^i^j^j^i[Vn+l,nGn.ffin,n+l + Vn+l,n+2G^j^^2ji+2pn+2,n+l] 

Now we require the terms G? „, G^^^ G^^^ „+2> G'n+2,2n+i= G'n+2,n+2- 
These can be found by using the above reduced Green's function technique, 
in addition to exploiting the recursive relation for the determinant of a matrix 
consisting only of diagonal and sub-diagonal entries. Similar calculational 
methods were employed by Evenson and Karplus [56]. 

For simplicity, we assume that the couplings to the two metals at the 
ends are the same, and use the wide-band approximation, in which A ~ 
7r^^t'^5(e — e/c) is taken as constant. 

If dn represents the determinant of a n x n matrix with diagonal entries 
set to e — ej and subdiagonal entries set to some then it is possible to 
express: 



i-v) 



n-1 



dn - dn[^] 



(8) 
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dn-l — dn-2[ 



dn — dn-l\ 



(9) 



6 



(10) 



^n+2,2n+l ~ ^l,n (H) 



Invoking the wide band approximation the above expressions reduce to 
the following form: 



^°""4 + kA ^^^^ 



Thus, the problem has been reduced to calculating the determinants dn, 
which is given in the appendix. 

The conductance can be obtained from Caroli's formula and integrated to 
obtained the net current. The net current thus obtained has a Qi, dependency 
which has to be eliminated by performing a thermal averaging. The final 
result obtained after thermal averaging gives the net total current. As noted 
in the introduction, in electrochemical systems there are two potential differ- 
ences to consider: the bias V between the two electrodes, and the electrode 
potential, which shifts the levels in the solution. We use the convention that 
the potential of the right electrode R is kept constant, and set its Fermi level 
to zero. The levels on the wire shift with the electrode potential; this assumes 
that the conductivity of the solution is higher than that of the wire. Other 
scenarios can be calculated by the same formalism. With this convention, we 
write the total current in the form: 

I{q) = J Tr[G'TiG"r«]{/(e + eoV) - f{e)}de (15) 
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where /(e) denotes the Fermi-Dirac distribution. As has been pointed out 
several times (see e.g. [54]), in the case of a classical solvent it is sufficient to 
consider a single effective solvent coordinate q. Effectively, this means that 
in the Hamiltonian, we make the following substitutions: 

^^fkJ^al Xq'^, ^hw^g^q^ -2Xq (16) 

The average over the solvent configurations can then be written as: 

7 = 1 y dge-^^(«)/(g) j dqe-^^^'^ (17) 

where the energy, as a function of the solvent coordinate g, is: 

E{q) ^Xq' + J2f 6(clc,) (18) 

Thus E{q) is obtained by performing a partial trace over the fermionic 
part of the total Hamiltonian. The quantity (c|cj), as viewed by Wingreen 
ct al. [18], is the lesser component of the Keldysh Green's function, G^^. 
If and Fji are the imaginary parts of the self-energy arising from the 
interaction with the left and the right reservoirs (which in the view of wide- 
band approximations is A), then for the present case 

G<^if{e + eoV)[GTLG% + tf{e)[G'rRG% (19) 

At this point, a few comments on the appearance for G^ are needed. 
It is well known that in equilibrium the lesser Green's function takes the 
form of a product of spectral function times the occupation function. (G^ = 
ia{e)f{e)). That is in case of zero bias V when both the reservoirs have the 
same potential, /(e -|- CqV) — /(e), then 

G< = tf{e)[G'-{rL + r«)G1 = z/(e)2A[G^G'^] (20) 

Now G^(e) — G"(e) = a(e) where a(e) is the spectral function, and the 
imaginary part of self-energy can be written as ( 1/G" — 1/ G*') . In our notation, 
2A = Al + = (1/G» - 1/G''). Also the diagonal part of G^(e) - G"(e) 
is proportional to the density of states, p{e). In equilibrium, we recover the 
result 

G< = z/(e)a(e) = /(e)(G'--G«). (21) 
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Substituting the above result in the expression (fTSi) for E{q), it is seen 
that at equihbrium E{q) = + J e/(e)ImTrG'(e)(ie, wherein the second 
term in the energy expression is widely employed in a variety of contexts in 
physics. Thus the expression for E{q) for the non-equilibrium case can be 
written compactly as 



E{q) = \q^ + J e/(e + eoV)TT[G'- rLG^]de + J e/(e)Tr[G'T^G"]rfe (22) 

The trace in the above equation runs over the (2n + 1) sites numbered 
by the index i. Now for better exposition of the computation involved in 
calculating the E{q), we consider a generic term which has the form shown 
below: 



E.(g) = A r eGi;iGt,rf6 + A [ eG',^^^,G^^^,^,de (23) 

J —OO J —CXI 

where we have replaced the Fermi-Dirac distribution by step function, and 
taken the Fermi level as zero. 

We consider the case for the three possible locations of i, {i < n,i > 
n + 2,i = n+1). As before the general idea behind the approach to get 
the matrix elements of G is to resort to a Dyson expansion. The closed form 
equation so obtained has to be solved to get the relevant terms. The choice 
of V is same as used before. 

Gl,.j = G;^ j + G;^ „'U„^„+lG„+l^i (24) 

G„+l,i = G°_,_;^^j + G°_^]^^„_^;^'U„+i_„+2G„+2,j + G^n+l,n+l^n+l,nGn,i (25) 
G„+2,i = G°_,_2,i + Gn+2,n+2^n+2,n+lGn+l,i (26) 
Gn,i = G°_j + G^„^Vn^n+lGn+l,i (27) 

Case I: i < n 

pO pO pO 

(28) 

dn + «Ad„_i 

g:;.^'-"'';'^--'^"'-^) (30) 
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Case II: i > n + 2 

no r"0 po 

'^l,n^n,n+l'^n+l,n+l^n+l,n+2»^n+2,i 



Gi, 



(31) 

It is crucial to note at this stage that because of the form of perturbation 
selected the unperturbed G° has a symmetric structure with respect to the 
first and third block matrix and hence Gn+2,i above is same as G° ^ in 

case I. 

Case III: i = n + 1 

^l,n ^n,n+l^n+l,n+l 



Gl,n 



1 ~ [G°+i,n+it'n+l,n+2G°+2,n+2'^n+2,n+l + Gn+l,n+l'^n+l,nGn,n'^n,n+l] 

(32) 

Even though i^n,n+i = '^n+i,n+2 = V, we have maintained the subscript 
indices for ease of checking the final expressions. 



3 Results and discussions 

The principle new feature of our work is the dynamic interaction of the 
redox system with the adjoining species, which fluctuates with the solvent 
coordinate q. The main effects can be demonstrated with a chain of three 
atoms, and we limit our numerical calculation to this case. 

Even though we have restricted our treatment to the symmetric case, in 
which the coupling A to the two leads and the interatomic couplings v are the 
same on both sides, the system contains a fair amount of parameters. In the 
following model calculations, we have set A = 0.3 eV and the reorganization 
energy A = 0.3 eV unless otherwise mentioned, and for the other parameters 
we have chosen values appropriate to demonstrate special effects. 

The case of a single intervening redox center is well examined. The new 
feature of the three-atom chain is the interaction between the levels ei of the 
two side atoms and the redox center. Before considering this in detail, it is 
instructive to investigate the reference case in which this effect is weak; in 
this limit, we should obtain similar results to the case of one atom. 
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The interaction of the side atoms with the center is weak, if e\ hes far 
from the Fermi level and the coupling v is small. In this limit, the redox 
center is at the equilibrium potential for = A. As expected, in this case 
the potential energy surfaces Ei^q) are similar to the one atom case [57]. 
At zero bias, they have the same form as for a normal, outer sphere redox 
reaction- see Fig. [TJ In the left well, the occupancy (n^,) is zero, in the right 
well unity. Application of a bias produces a region with (ra^) ^1/2, which 
extends the barrier in the center. For the three atom case, the tunneling rate, 
as a function of the electronic energy e and the solvent coordinate g, is given 



As long as ei lies so high that it plays no role, this rate has a maximum 
where e — + 2Ag = 0. Tunneling occurs only between the two Fermi levels, 
in the range —CqV < e < 0. Inspection shows, that the maximum of t(e, q) 
is obtained in the region where (n^) ~ 1/2, which therefore gives the main 
contribution to the current. 

Continuing with the case of weak coupling, the current-potential curves 
are symmetric at the equilibrium condition = A. Shifting by application 
of an overpotential leads to asymmetrical curves with rectifying properties. 
Figure [2] shows the case in which the redox level has been lowered, so that 
the most favorable energy range now lies below the Fermi level of the right 
electrode. Therefore the current is higher at positive bias, where this energy 
range lies between the two Fermi levels. 

Really new features occur when the redox level interacts noticeably with 
the levels ei. In general, three interacting atomic levels combine to form three 
molecular orbitals. Since the redox level changes its energy with the solvent 
fluctuations, so do the resulting molecular orbitals. So, for some range of q the 
redox level will be far from ei and the interaction will be almost negligible, 
in another range it will lie close in energy, so that one observes the typical 
splitting of the levels. 

The transition probability t{e, q) has local maxima, whenever e is near one 
of the molecular orbitals. Since is depends only on the combination e,. — 2Ag, 
it is sufficient to investigate the dependence on e,. for g = 0, as is done 
in Fig. [31 The interesting region lies where ~ e^. When both are equal, 
there are three distinct maxima at the three molecular orbitals. When they 
are separated, only two maxima occur at the position of the atomic orbitals. 



by: 




(e - e, + 2Ag)2 [(e - + 2Ag)(e - d + 2 A) - 2v'^f 



(33) 
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since the splitting is too small to show up - it is hidden beneath the maxima. 

These oscillations in the transition probability give rise to interesting 
current-potential curves exhibiting several steps and even regions with a neg- 
ative differential resistance (see Fig. Hj), effects which do not occur with a 
single electronic intermediate state. The exact shape of these characteristics 
is determined by an interplay of three effects: The change of the potential- 
energy curves with the bias, the dynamic changes in the energy of the molec- 
ular orbitals as the solvent coordinate q fluctuates, and the resulting oscil- 
lations in the transition rate. These highly nonlinear effects are more pro- 
nounced when the coupling A to the two leads is weaker. Figure [5] shows two 
examples where the system parameters have been chosen such that the curves 
either exhibit nice plateaus or a pronounced negative differential resistance. 

As pointed out in the introduction, in electrochemical systems two volt- 
ages can be controlled independently, the bias and the potential between the 
solution and one electrode. This makes it possible to perform spectroscopy 
of the electronic states in the wire, which experience the potential of the 
solution. We introduce the overpotential rj of the redox couple with respect 
to the right electrode through = A — eo?7, and let the intermediate state 
shift in the same way: ei = — cq?]. In a real system, because of the finite 
conductivity of the solution, 77 may be only a fraction of the externally ap- 
plied potential, but this would require only a trivial modification. If we keep 
the bias constant at a comparatively small value and scan the overpotential 
?7, we obtain a peak in the current every time an electronic state lies within 
the tunneling range of energy between the two Fermi levels. A few examples 
are seen in Fig. [61 The redox level always gives a peak near 77 = 0, and for 
the parameters chosen we see a second peak near e?. When these two en- 
ergies lie close, one peak may appear as a shoulder. Note that the curves 
for = ±0.5 eV in the figure are not quite symmetric, because the bias 
breaks the symmetry. In theory, we could expect to see up to three peaks in 
these curves corresponding to the three molecular orbitals formed, but these 
only occur for a very strong coupling v and small energies of reorganisation. 
Otherwise the splitting induced by v is hidden under the peak for ei. 

4 Conclusions 

In this work we have presented a model for the conductivity of a molecular 
wire containing a redox system, and embedded in an electrochemical envi- 
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ronment. We considered the interaction of the redox system with a classical 
solvent, whose state was represented by a solvent coordinate q in the spirit 
of the Marcus theory. Using the wide-band approximation, we were able to 
derive an exact expression for the quantum conductance of a chain of arbi- 
trary length. The thermal average over the solvent configurations had to be 
performed numerically. 

Explicit calculations have been performed for a chain of three atoms. 
When the electronic levels of the neighboring atoms interact weakly with 
the redox couple - because their energies are very different or the coupling is 
weak - the wire behaves much like a single intervening atoms. Interesting new 
features arise when the redox couple interacts strongly with the neighboring 
levels. Since the redox level fluctuates with the solvent, this interaction is 
dynamic and changes with the solvent configuration. This gives rise to inter- 
esting relations between the current through the wire and the applied bias. In 
particular, rectification, extended current plateaus, and negative differential 
resistances can be observed. 

In an electrochemical environment two potential differences, the bias and 
the electrochemical potential of the wire, can be varied independently. This 
makes it possible to perform electronic spectroscopy at constant bias by 
changing the electrochemical potential. Intermediate states show up as char- 
acteristic current peaks. 

Our treatment has been limited to a redox couple interacting with a 
classical solvent. An extension to the case where quantum modes couple to 
the electron transfer should be possible, using Green's function techniques 
that have been apphed to the case of a single atom [40] . This could give rise 
to additional structure in current-potential curves. 
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A Calculation of the determinants dn 

The method for obtaining the determinants dn is discussed in this section. 
Let the determinant of a n x n tri-diagonal matrix with diagonal entries e 
and sub-diagonal entries v be denoted as d„. Then the following recursion 
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holds good: 



V 



M 



V 

e 



\ 



v: 



V 



V 



(34) 



dn = edn-l " (35) 

The above recursive relation can be written compactly in a matrix form: 



dn 

dfi-i 



e —V 
1 



dn-2 



(36) 



Denoting the recursion matrix as R, it's eigenvalues Ai and A2 and the matrix 
S which diagonalises R can be found easily. Iterating the recursion relation 
while noting that do — 1 and di — e, gives the following result: 



n-l 



Ai 
A, 



n-l 



di 
do 



(37) 



The final required expression for dn obtained from the above equation is 
shown below: 



dn = 7 7 With Ai,2 = 

Ai — A2 z 



(38) 



The methods discussed above can be use d to get the inverse elements of 
the matrix: 



dn-i 



—vd. 



{-vY'^dr, 



dn 



'12 



dn 



dr, 



(39) 
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Figures 




Figure 1: Potential energy surfaces in the case of weak coupling and at the 
equilibrium potential for the redox system; system parameters: ei = 0.8 eV, 
V = 0.01 eV, = A = 0.3 eV. The insert shows the occupation (ra^) of the 
redox center. 
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Figure 2: Current-potential curves in the case of weak coupling; system pa- 
rameters: ei = 1.2 eV, V = 0.01 eV, 




Figure 3: Transition probability as a function of the energy e of the tunneling 
electron for various values of e^; the solvent coordinate q was set to zero. 
System parameters: : ei = 0.2 eV, f = 0.1 eV. 
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Figure 4: Current-potential curves for various values of e^. System parame- 
ters: : ei = 0.2 eV, v = 0.1 eV. 



Figure 5: Current-potential curves for small coupling to the leads. System 
parameters: A = 0.1 eV, = 0.5 eV, v = 0.1 eV, ei = —0.2 eV (left curve) 
and ei = 0.3 eV (right curve). 
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Figure 6: Current at constant bias V = 0.1 V as a function of the overpo- 
tential 1] for various values of e?, the value of ei for vanishing overpotential; 
V = 0.1 eV. 
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